
****************************************************************************************************
*Table A.6: Effects on quantity disbursed, quantity received, and leakage
*This table reports estimated treatment effects on the quantity of commodities disbursed by the government, received by recipients, and the difference in endline one (January - March)
****************************************************************************************************

use "$SurveyDataDir/JH_ePOS_HH_DataforAnalysis.dta",clear

* Keep if surveyed or classified as ghost
keep if surveyed == 1 | ghost_final == 1                 


svyset [pw = pweight]


count if ghost_final == 1
scalar ghosts = r(N)
scalar obs = 3960   

*Get relative weights of AAY and PH rationcard holders by whether RC is in an urban area
sum pweight if rationcardtype == "AAY" & isurban == 0
scalar AAY_weight0 = r(sum)

sum pweight if rationcardtype == "PH" & isurban == 0
scalar PH_weight0 = r(sum)

sum pweight if rationcardtype == "AAY" & isurban == 1
scalar AAY_weight1 = r(sum)

sum pweight if rationcardtype == "PH" & isurban == 1
scalar PH_weight1 = r(sum)


***************************************************************
*store mean entitlement in place holder regression results

loc ration "rice wheat sugar salt kero"

foreach rat of local ration {

*Estimation
preserve

keep value_`rat'_mar17 value_`rat'_feb17 value_`rat'_jan17 treatment strata pweight block_code uid rationcardtype stat_ent_quantity_`rat'* rationcardtype isurban
gen value_`rat'3 = value_`rat'_mar17
gen value_`rat'2 = value_`rat'_feb17
gen value_`rat'1 = value_`rat'_jan17
drop value_`rat'_mar17 value_`rat'_feb17 value_`rat'_jan17

reshape long value_`rat', i(uid) j(month)

count
scalar obs_all = 3*3960
*Statutory NIC entitlement

qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "AAY" & isurban == 0
qui estat sd
matrix MeanAAY0s = r(mean)

qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "PH" & isurban == 0
qui estat sd
matrix MeanPH0s = r(mean)

qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "AAY" & isurban == 1
qui estat sd
matrix MeanAAY1s = r(mean)

qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "PH" & isurban == 1
qui estat sd
matrix MeanPH1s = r(mean)

di in red "`rat'"
matrix list MeanAAY0s
matrix list MeanAAY1s
matrix list MeanPH0s
matrix list MeanPH1s


scalar avg_stat_ent_quantity_`rat' = (AAY_weight0*`=MeanAAY0s[1,1]' + PH_weight0*`=MeanPH0s[1,1]' + AAY_weight1*`=MeanAAY1s[1,1]' + PH_weight1*`=MeanPH1s[1,1]')/(AAY_weight0 + PH_weight0 + AAY_weight1 + PH_weight1)

qui svy: mean value_`rat' if treatment == 0 
qui estat sd
matrix Mean0 = r(mean)

eststo `rat'_v1BL: xi: reg value_`rat' treatment i.strata [pw = pweight], cluster(block_code)
estadd scalar stat_ent_mean = avg_stat_ent_quantity_`rat'

restore


}

di in red "EL1"


loc tabname "${OutputDir}/TableA_5.tex" 
cd "${adoDir}"
				
MultiPartTabStart, ///
			ncol(5) tabname(`tabname') ///
			colnames("Rice" "Wheat" "Sugar" "Salt" "Kerosene") width("\hsize")			

**** Mean Entitlement 

MultiPartTabPanelStartEntitle, ///
			ncol(5) tabname(`tabname') ///
			panelstring("")
			
			estimates restore *rice_v1BL
			estimates restore *wheat_v1BL
			estimates restore *sugar_v1BL
			estimates restore *salt_v1BL
			estimates restore *kero_v1BL
			
MultiPartTabPanelEndEntitle, ///
				ncol(6) tabname(`tabname') ///
				models("*rice_v1BL *wheat_v1BL *sugar_v1BL *salt_v1BL *kero_v1BL") ///
				drop(treatment _Ist* _cons) ///
				cells(b) ///
				stats(stat_ent_mean, label("\textit{Mean entitlement}") fmt(0))
	


	
	
	
	
**************************************************************
* 1)disbursement
**************************************************************

********************************************
* Use EL1 disbursement data
********************************************
use "${AdminDataDir}/allocation_FPSlevel_jan17_mar17.dta",clear

**Average entitlement quantity (listed in summary stats portion of table)

local ration "rice wheat salt sugar kero"


*******
* Pooled
*******
local tflist ""
local tflist_noBL ""
local tflist_trend ""
local tflist_trend_noBL "" 

foreach rat of local ration{
scalar obs=8924*3


gen bl_var = dis_qty_perRC_`rat'_y0
egen bl_var_mean = mean(bl_var)
gen bl_var_mi = missing(bl_var)
replace bl_var = bl_var_mean if bl_var_mi == 1



preserve
keep dis_qty_perRC_`rat'_*17 ent_qty_perRC_`rat'_*17 rc_count_*17 treatment strata block_code bl_var bl_var_mi admin_dealer_id
rename dis_qty_perRC_`rat'_jan17 dis_qty_perRC1
rename dis_qty_perRC_`rat'_feb17 dis_qty_perRC2
rename dis_qty_perRC_`rat'_mar17 dis_qty_perRC3
rename ent_qty_perRC_`rat'_jan17 ent_qty_perRC1
rename ent_qty_perRC_`rat'_feb17 ent_qty_perRC2
rename ent_qty_perRC_`rat'_mar17 ent_qty_perRC3
rename rc_count_jan17 rc_count1
rename rc_count_feb17 rc_count2
rename rc_count_mar17 rc_count3

reshape long dis_qty_perRC ent_qty_perRC rc_count, i(admin_dealer_id) j(month)

gen treatmentXmonth=treatment*month

qui: summarize dis_qty_perRC if treatment == 0 [aw=rc_count]
matrix Mean0 = r(mean)

qui: mean ent_qty_perRC if treatment == 0 
qui estat sd
matrix Mean0e = r(mean)


eststo pool_`rat'_q1BL: xi: reg dis_qty_perRC treatment bl_var bl_var_mi i.strata[aweight=rc_count] , cluster(block_code) 
estadd scalar control_mean = `=Mean0[1,1]'
estadd scalar stat_ent_mean = `=Mean0e[1,1]'
estadd scalar percent_obs=100*e(N)/obs


tempfile pool_`rat'_q1BL
parmest, label saving("`pool_`rat'_q1BL'")
local tflist "`tflist' `pool_`rat'_q1BL'"


restore

drop bl_*
}

*Pooled with baseline

loc ration "rice wheat sugar salt kero"
dsconcat `tflist',subset(parm p in 1) 
qqvalue p, method(yekutieli) qvalue(qval)
mkmat qval

*MHT-adding FDR adjusted p values to stored estimates
forval i = 1/5 {
local rat: word `i' of `ration'
estimates restore pool_`rat'_q1BL


mat bet=e(b)                   

mat bet[1,1]=qval[`i',1]             
mat list bet      
estadd matrix qvalue = bet               
                                            
mat list e(qvalue)
}



	
loc nc=5
			
MultiPartTabPanelStart, ///
			ncol(`nc') tabname(`tabname') ///
			panelstring("Panel A: Quantity disbursed")
			
MultiPartTabPanelEnd, ///
				ncol(`nc') tabname(`tabname') ///
				models("*rice_q1BL *wheat_q1BL *sugar_q1BL *salt_q1BL *kero_q1BL") ///
				drop(_cons _Istrata* bl_var_mi bl_var) ///
				cells(b(star fmt(3) ) se(par(( )) fmt(3) ) qvalue(par([ ])  keep(treatment) fmt(2)) )  ///
				varlabels(treatment "Treatment", elist(treatment \addlinespace )  ) ///
				stats(control_mean N, label("Control mean" "Observations") fmt(%12.2g %15.0fc)) starlevels( * 0.10 ** 0.05 *** 0.01)


**************************************************************
* 1) Receipt
**************************************************************

use "${SurveyDataDir}/JH_ePOS_HH_DataforAnalysis.dta",clear

* Keep if surveyed or classified as ghost
keep if ss_code == "SS01" | ghost_final == 1

svyset [pw = pweight]


count if ghost_final == 1
scalar ghosts = r(N)
scalar obs = 3960   

*Get relative weights of AAY and PH rationcard holders by whether RC is in an urban area
sum pweight if rationcardtype == "AAY" & isurban == 0
scalar AAY_weight0 = r(sum)

sum pweight if rationcardtype == "PH" & isurban == 0
scalar PH_weight0 = r(sum)

sum pweight if rationcardtype == "AAY" & isurban == 1
scalar AAY_weight1 = r(sum)

sum pweight if rationcardtype == "PH" & isurban == 1
scalar PH_weight1 = r(sum)



local tflist ""
local tflist_noBL ""
local tflist_trend ""
local tflist_trend_noBL ""
foreach rat of local ration {


gen bl_var = `rat'_weigh_y0
egen bl_var_mean = mean(`rat'_weigh_y0)
gen bl_var_mi = missing(`rat'_weigh_y0)
replace bl_var = bl_var_mean if bl_var_mi == 1



preserve
keep uid b16_`rat'_weigh_jan17 b16_`rat'_weigh_feb17 b16_`rat'_weigh_mar17 ent_quantity_`rat'_* stat_ent_quantity_`rat' bl_var bl_var_mi treatment strata pweight block_code rationcardtype isurban
gen b16_`rat'_weigh3 = b16_`rat'_weigh_mar17
gen b16_`rat'_weigh2 = b16_`rat'_weigh_feb17
gen b16_`rat'_weigh1 = b16_`rat'_weigh_jan17
gen ent_quantity_`rat'3 = ent_quantity_`rat'_mar17
gen ent_quantity_`rat'2 = ent_quantity_`rat'_feb17
gen ent_quantity_`rat'1 = ent_quantity_`rat'_jan17
drop b16_`rat'_weigh_mar17 b16_`rat'_weigh_feb17 b16_`rat'_weigh_jan17 ent_quantity_`rat'_mar17 ent_quantity_`rat'_feb17 ent_quantity_`rat'_jan17
reshape long b16_`rat'_weigh ent_quantity_`rat', i(uid) j(month)

count
scalar obs_all = 3*3960

gen treatmentXmonth=treatment*month


qui svy: mean b16_`rat'_weigh if treatment == 0 
qui estat sd
matrix Mean0 = r(mean)

* Food department monthly entitlement
qui svy: mean ent_quantity_`rat' if rationcardtype == "AAY" & isurban == 0
qui estat sd
matrix MeanAAY0 = r(mean)

qui svy: mean ent_quantity_`rat' if rationcardtype == "PH" & isurban == 0
qui estat sd
matrix MeanPH0 = r(mean)

qui svy: mean ent_quantity_`rat' if rationcardtype == "AAY" & isurban == 1
qui estat sd
matrix MeanAAY1 = r(mean)

qui svy: mean ent_quantity_`rat' if rationcardtype == "PH" & isurban == 1
qui estat sd
matrix MeanPH1 = r(mean)

*Statutory NIC entitlement
qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "AAY" & isurban == 0
qui estat sd
matrix MeanAAY0s = r(mean)

qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "PH" & isurban == 0
qui estat sd
matrix MeanPH0s = r(mean)

qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "AAY" & isurban == 1
qui estat sd
matrix MeanAAY1s = r(mean)

qui svy: mean stat_ent_quantity_`rat' if rationcardtype == "PH" & isurban == 1
qui estat sd
matrix MeanPH1s = r(mean)

eststo `rat'poolBL1: xi: reg b16_`rat'_weigh treatment bl_var bl_var_mi i.strata [pw = pweight], cluster(block_code)
estadd scalar control_mean = `=Mean0[1,1]'
estadd scalar percent_obs 100*e(N)/obs_all
estadd scalar ent_mean = (AAY_weight0*`=MeanAAY0[1,1]' + PH_weight0*`=MeanPH0[1,1]' + AAY_weight1*`=MeanAAY1[1,1]' + PH_weight1*`=MeanPH1[1,1]')/(AAY_weight0 + PH_weight0 + AAY_weight1 + PH_weight1)
estadd scalar stat_ent_mean = (AAY_weight0*`=MeanAAY0s[1,1]' + PH_weight0*`=MeanPH0s[1,1]' + AAY_weight1*`=MeanAAY1s[1,1]' + PH_weight1*`=MeanPH1s[1,1]')/(AAY_weight0 + PH_weight0 + AAY_weight1 + PH_weight1)

tempfile `rat'poolBL1
parmest, label saving("``rat'poolBL1'")
local tflist "`tflist' ``rat'poolBL1'"



restore

drop bl_var*

}

dsconcat `tflist', subset(parm p in 1) dsn(ds_ration)  
qqvalue p, method(yekutieli) qvalue(qval)
mkmat qval

 
forval i = 1/5 {
local rat: word `i' of `ration'
estimates restore `rat'poolBL1

mat bet=e(b)                   

mat bet[1,1]=qval[`i',1]             
mat list bet      
estadd matrix qvalue = bet               
                                            
mat list e(qvalue) 
}



MultiPartTabPanelStart, ///
			ncol(`nc') tabname(`tabname') ///
			panelstring("Panel B: Quantity received")
			
MultiPartTabPanelEnd, ///
				ncol(`nc') tabname(`tabname') ///
				models("ricepoolBL1 wheatpoolBL1 sugarpoolBL1 saltpoolBL1 keropoolBL1") ///
				drop(_cons _Istrata* bl_var_mi bl_var) ///
				cells(b(star fmt(%12.2g) ) se(par(( )) fmt(%12.2g) ) qvalue(par([ ])  keep(treatment) fmt(2)) )  ///
				varlabels(treatment "Treatment", elist(treatment \addlinespace )  ) ///
				stats(control_mean N, label("Control mean" "Observations") fmt(%12.2g %15.0fc)) starlevels( * 0.10 ** 0.05 *** 0.01)
				

				
				
**************************************************************
* 3) Leakage
**************************************************************
			
		
********************************************
* Use EL1 disbursement data
********************************************
use "${AdminDataDir}/allocation_blocklevel_jan17_mar17.dta",clear

preserve
use "${SurveyDataDir}/JH_ePOS_HH_DataforAnalysis.dta",clear
*Keep only variables needed in analysis
keep if ss_code == "SS01" | ghost_final == 1 

keep ss_code ghost_final district_code block_code fps_code pweight uid treatment isurban strata b16*17 value* *weigh*_y0
tempfile HHdata
save `HHdata'
restore

*merge with household data
merge 1:m district_code block_code isurban strata treatment using `HHdata',update  

di in red "`ration'"

svyset block_code [pw=pweight]


loc rat "kero"
local tflist ""

foreach rat of local ration{
preserve
	*BL vars
	gen bl_var_d = dis_qty_perRC_`rat'_y0
	egen bl_var_mean_d = mean(bl_var_d)
	gen bl_var_mi_d = missing(bl_var_d)
	replace bl_var_d = bl_var_mean_d if bl_var_mi_d == 1

	gen bl_var_r = `rat'_weigh_y0
	egen bl_var_mean_r = mean(bl_var_r)
	gen bl_var_mi_r = missing(bl_var_r)
	replace bl_var_r = bl_var_mean_r if bl_var_mi_r == 1
	
	*Reshape
	keep uid b16_`rat'_weigh_*17 dis_qty_perRC_`rat'_*17 ent_qty_perRC_`rat'_*17 rc_count_*17 bl_var_r bl_var_mi_r bl_var_d bl_var_mi_d treatment strata pweight block_code
	rename b16_`rat'_weigh_mar17 b16_`rat'_weigh3
	rename b16_`rat'_weigh_feb17 b16_`rat'_weigh2
	rename b16_`rat'_weigh_jan17 b16_`rat'_weigh1
	rename dis_qty_perRC_`rat'_jan17 dis_qty_perRC1
	rename dis_qty_perRC_`rat'_feb17 dis_qty_perRC2
	rename dis_qty_perRC_`rat'_mar17 dis_qty_perRC3
	rename ent_qty_perRC_`rat'_jan17 ent_qty_perRC1
	rename ent_qty_perRC_`rat'_feb17 ent_qty_perRC2
	rename ent_qty_perRC_`rat'_mar17 ent_qty_perRC3
	rename rc_count_jan17 rc_count1
	rename rc_count_feb17 rc_count2
	rename rc_count_mar17 rc_count3
	
	reshape long b16_`rat'_weigh dis_qty_perRC ent_qty_perRC rc_count, i(uid) j(month)
	
	count
	scalar obs_all = 3*3960
	
	gen treatmentXmonth=treatment*month
	
	*Disbursement per RC
	qui svy: mean dis_qty_perRC if treatment == 0 
	qui estat sd
	matrix Mean0d = r(mean)
	
	qui mean ent_qty_perRC if treatment == 0
	qui estat sd
	matrix Mean0e = r(mean)
	
	
	count if (!missing(dis_qty_perRC) & !missing(b16_`rat'_weigh))   
	scalar used_obs=r(N)  
	svyset block_code [pw=pweight]

	
	*Disbursement 
	eststo `rat'_dis_pooled: xi: svy:reg dis_qty_perRC treatment bl_var_d bl_var_mi_d i.strata 
	estadd scalar control_mean = `=Mean0d[1,1]'
	estadd scalar stat_ent_mean = `=Mean0e[1,1]'
	
	*Receipt
	qui svy: mean b16_`rat'_weigh if treatment == 0 
	qui estat sd
	matrix Mean0r = r(mean)
	
	eststo `rat'_rep_pooled: xi: svy:reg b16_`rat'_weigh treatment bl_var_r bl_var_mi_r i.strata 
	estadd scalar control_mean = `=Mean0r[1,1]'
	estadd scalar stat_ent_mean = `=Mean0e[1,1]'
	
	
	*Jointly estimate leakage
	suest `rat'_dis_pooled `rat'_rep_pooled
	scalar r2_adj = e(r2_a)
	
	eststo pooled_`rat': lincomest [`rat'_dis_pooled]treatment - [`rat'_rep_pooled]treatment
	estadd scalar control_mean = `=Mean0d[1,1]' - `=Mean0r[1,1]'
	estadd scalar used_obs=used_obs
    estadd scalar percent_obs 100*used_obs/obs_all
	estadd scalar r2_adj = r2_adj
	estadd scalar stat_ent_mean = `=Mean0e[1,1]'
	
	tempfile pooled_`rat'
	parmest, label saving("`pooled_`rat''")
	local tflist "`tflist' `pooled_`rat''"
	
restore
}


preserve
dsconcat `tflist', subset(parm p in 1) dsn(ds_ration)  
qqvalue p, method(yekutieli) qvalue(qval)
mkmat qval


loc ration_all "rice wheat sugar salt kero"
forval i = 1/5{
local estname: word `i' of `tflist'
local rat: word `i' of `ration_all'

estimates restore pooled_`rat'
mat bet=e(b)                   

mat bet[1,1]=qval[`i',1]             
mat list bet      
estadd matrix qvalue = bet               
                                            
mat list e(qvalue) 
}
			
				
#delimit cr
				
MultiPartTabPanelStart, ///
			ncol(`nc') tabname(`tabname') ///
			panelstring("Panel C: Leakage")
			
MultiPartTabPanelEnd, ///
				ncol(`nc') tabname(`tabname') ///
				models("pooled_rice pooled_wheat pooled_sugar pooled_salt pooled_kero") ///
				cells(b(star fmt(%12.2g) ) se(par(( )) fmt(%12.2g) ) qvalue(par([ ]) fmt(2)) )  ///
				varlabels((1) "Treatment", elist((1) \addlinespace )  ) ///
				stats(control_mean used_obs, label("Control mean" "Observations") fmt(%12.2g %15.0fc)) starlevels( * 0.10 ** 0.05 *** 0.01)
				
			
MultiPartTabEnd, ///
			ncol(`nc') tabname(`tabname')
